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(57) Abstract: The discrete Fourier transform (DFT) 
is computed in a plurality of parallel processors. A 
DFT of length N is divided into r partial DFTs of length 
(N/r), in which the r partial DFTs are calculated in 
separate parallel processors and then combined in a 
combination phase to form a complete DFT of length 
(N). The r partial FFTs are able to be computed in 
parallel multiprocessors by defining the mathematical 
model of the combination phase in such manner so as to 
allow the r parallel processors to operate independently 
and simultaneously. A second embodiment presents 
a radix-r fast Fourier algorithm that reduces the 
computational effort as measured by the number of 
multiplications and permits the N/r parallel processors 
to operate simultaneously and with a single instruction 
sequence. 
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PARALLEL MULTIPROCESSING FOR THE FAST FOURIER TRANSFORM WITH 
PIPELINE ARCHITECTURE 

Field of Invention 

The present invention relates to the field of Fast Fourier Transform analysis. In particular, 
the present invention relates to a parallel processing architecture adapted for use in a 
pipelined Fast Fourier Transform method and apparatus. 

Background of the invention 

Physical parameters such as light, sound, temperature, velocity and the like are converted 
to electrical signals by sensors. An electrical signal may be represented in the time domain 
as a variable that changes with time. Alternatively, a signal may be represented in the 
frequency domain as energy at specific frequencies. In the time domain, a sampled data 
digital signal is a series of data points corresponding to the original physical parameter. In 
the frequency domain, a sampled data digital signal is represented in the form of a plurality 
of discrete frequency components such as sine waves. A sampled data signal is 
transformed from the time domain to the frequency domain by the use of the Discrete 
Fourier Transform (DFT). Conversely, a sampled data signal is transformed back from the 
frequency domain into the time domain by the use of the Inverse Discrete Fourier 
Transform (IDFT). 

The Discrete Fourier Transform is a fundamental digital signal-processing transformation 
used in many applications. Frequency analysis provides spectral information about signals 
that are further examined or used in further processing. The DFT and IDFT permit a signal 
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to be processed in the frequency domain. For example, frequency domain processing 
allows for the efficient computation of the convolution integral useful in linear filtering 
and for signal correlation analysis. Since the direct computation of the DFT requires a 
large number of arithmetic operations, the direct computation of the DFT is typically not 
used in real time applications. 

Over the past few decades, a group of algorithms collectively known as Fast Fourier 
Transform (FFT) have found use in diverse applications, such as digital filtering, audio 
processing and spectral analysis for speech recognition. The FFT reduces the 
computational burden so that it may be used for real-time signal processing. In addition, 
the fields of applications for FFT analysis are continually expanding. 

Computational burden 

Computation burden is a measure of the number of calculations required by an algorithm. 
The DFT process starts with a number of input data points and computes a number of 
output data points. For example, an 8-point DFT may have an 8-point output. See figure 
1A. The DFT function is a sum of products, i.e., multiplications to form product terms 
followed by the addition of product terms to accumulate a sum of products (multiply- 
accumulate, or MAC operations). See equation (1) below. The direct computation of the 
DFT requires a large number of such multiply-accumulate mathematical operations, 
especially as the number of input points is made larger. Multiplications by the twiddle 
factors Wx, dominate the arithmetic workload. 
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To reduce the computational burden imposed by the computationally intensive DFT, 
previous researchers developed the Fast Fourier Transform (FFT) algorithms in which the 
number of required mathematical operations is reduced. In one class of FFT methods, the 
computational burden is reduced based on the divide-and-conquer approach. The principle 
of the divide-and-conquer approach method is that a large problem is divided into smaller 
sub-problems that are easier to solve. In the FFT case, the division into sub-problems 
means that the input data are divided in subsets for which the DFT is computed to form 
partial DFTs. Then the DFT of the initial data is reconstructed from the partial DFTs. See 
N. W. Cooley and J. W. Tukey, " An algorithm for machine calculation of complex 
Fourier series", Math.Comput., Vol. 19 pp. 297-301, Apr 1965. There are two approaches 
to dividing (also called decimating) the larger calculation task into smaller calculation sub- 
tasks: decimation in frequency (DIF) and decimation in time (DIT). 

Butterfly implementation of the DFT 

For example, an 8-point DFT can be divided into four 2-point partial DFTs. The basic 2- 
point partial DFT is calculated in a computational element called a radix-2 DIT butterfly 
(or butterfly-computing element) as represented in figure 2A1. Similarly to the DIT 
butterfly-computing element, figure 2A2 shows the function of a radix-2 DIF butterfly. A 
radix-2 butterfly has 2 inputs and 2 outputs, and computes a 2-point DFT. Figure 2B shows 
an FFT using 12 radix-2 butterflies to compute an 8-point DFT. Butterfly-computing 
elements are arranged in stages. There are three stages 1302, 1304 and 1306 of butterfly 
calculation. Data, x n is fed to the input of the butterfly-computing elements in the first 
stage 1302. After the first stage 1302 of butterfly-computation is complete, the result is fed 
to the in input of the next stage(s) of butterfly-computing element(s) and so on. 
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In particular, four radix-2 butterflies operate in parallel in the first stage 1302 to compute 8 
partial DFTs. The 8 outputs of the first stage 1302 are combined in 2 additional stages 
1304, 1306 to form a complete 8-point DFT output, X n . Specifically, the second stage 
1304 of 4 radix-2 butterflies and the third stage 1306 of 4 radix-2 butterflies comprise a 
two stage combination phase in which 8 radix-2 butterflies responsive to 8 partial DFTs 
form the final 8-point DFT function, X n . 

Figure 2C shows an FFT using 32 radix-2 butterflies to compute a 16-point DFT. There are 
4 stages of butterfly calculation. Eight radix-2 butterflies operate in parallel in the first 
stage 1402 where 2-point partial DFTs are calculated. The outputs of the first stage are 
combined in 3 additional combination stages 1403, 1404 and 1406 to form a complete 16- 
point DFT output. The output of the second stage 1403 of 8 radix-2 butterflies is coupled 
to a third stage 1404 of 8 radix-2 butterflies. The output of the third stage 1404 of 8 radix-2 
butterflies is coupled to a fourth stage 1406 of 8 radix-2 butterflies, the output of which the 
final 16-point DFT function. The combination phases 1403, 1404, 1406 comprise a 
combination phase in which 24 radix-2 butterflies responsive to 16 partial DFTs (from the 
first phase 1402) form the final 16 point DFT function, X n . 

Higher order butterflies may be used. See figure 2D, which uses 8 radix-4 butterflies in 2 
stages 1502, 1502 to compute a 16-point DFT. In general, a radix-r butterfly is a 
computing element that has r input points and calculates a partial DFT of r output points. 
In figure 2D, four radix-4 butterflies compute 16 partial DFTs in a first stage 1502. The 
combination phase 1504 comprises four radix-4 butterflies responsive to 16 partial DFTs 
(from the first phase 1502) to form the final 16 point DFT function, X n . 
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Communication burden 

A computational problem involving a large number of calculations may be performed one 
calculation at a time by using a single computing element. While such a solution uses a 
minimum of hardware, the time required to complete the calculation may be excessive. To 
speed up the calculation, a number of computing elements may be used in parallel to 
perform all or some of the calculations simultaneously. A massively parallel computation 
will tend to require an excessively large number of parallel-computing elements. Even so, 
parallel computation is limited by the communication burden. For example, a large number 
of data and constants may have to be retrieved from memory over a finite capacity data 
bus. In addition, intermediate results in one parallel-computing element may have to be 
communicated to another parallel-computing element. The communication burden of an 
algorithm is a measure of the amount of data that must be moved, and the number of 
calculations that must be performed in sequence (i.e., that cannot be performed in parallel). 

In particular, in a butterfly implementation of the DFT, some of the butterfly calculations 
cannot be performed simultaneously, i.e., in parallel. Subsequent stages of butterflies 
cannot begin calculations until earlier stages of butterflies have completed prior 
calculations. Also, the connections between butterflies in each stage to butterflies in the 
other stages impose a heavy communication burden between the butterfly computation 
stages. Thus, parallel implementations of the butterfly DFT are hampered by a heavy 
communication burden between butterflies. 



Summary of the invention 
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The heavy communication burden between the butterfly stages in the prior art results from 
structuring the butterfly implementation such that the first butterfly stage computes partial 
DFTs over the input data, and the latter butterfly stages combine the partial DFTs. In 
accordance with the present invention, partial DFTs are computed in a plurality of separate 
parallel processors and then combined in a single stage of combination for the parallel 
processing algorithm. Also, In accordance with the present invention for the multi-stage 
parallel processing algorithm, partial DFTs are computed in a plurality of separate parallel 
circuit boards that contain plurality of separate parallel chips, which contain plurality of 
separate parallel processors. The output data is obtained by combining firstly the outputs of 
the plurality of separate parallel processors, secondly by combining the outputs of the 
plurality of separate parallel chips and finally by combining the plurality of the outputs of 
the separate parallel circuit boards. 

The present architecture is a reorganization of the butterfly calculation of the DFT so as to 
reduce the communication burden between butterflies implemented in parallel computing 
elements. In particular, no communication is required between pluralities of separate 
parallel processors, chips or circuit boards. A combination stage of butterfly calculation is 
provided which combines the outputs of all the parallel processors (chips or circuit 
boards). 

In accordance with the present invention, the input data points of an N point DFT are 
divided into subsets of data. A plurality of processors, chips or circuit boards operate in 
parallel and independently of each other, and each being responsive to each respective 
subset of the input data points. The partial DFTs at the output of the plurality of parallel 
processors in each chip are then combined in a single combination phase to provide the 
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complete partial DFT for a single chip. The partial DFTs at the output of the plurality of 
parallel chips in each circuit board are then combined in a single combination phase to 
provide the complete partial DFT for a single board. The partial DFTs at the output of the 
plurality of parallel circuit boards are then combined in a single combination phase to 
provide the DFT solution. 

In the general case, each output data point of the DFT is a function of all of the input data 
points. However, by dividing the input data set into subsets, and operating each parallel 
processor independently of the other parallel processors in accordance with the present 
invention, the communication between the parallel processors is eliminated, thereby 
reducing the communication burden. Each or several parallel processors may then be 
implemented on a separate semiconductor chip or circuit board, without requiring any 
communication with any of the other parallel processors. 

In a second embodiment of the present invention, a radix-r butterfly implementation is 
provided in which the plurality of independent processors are operated in parallel using the 
same instructions and accessing the same necessary set of multiplier coefficients from 
memory at the same time. The resulting algorithm, in which a number of parallel 
processors operate simultaneously by a single instruction sequence, reduces both the 
computational burden and the communication burden. 

Brief description of the drawings 

Figure 1 A is a block diagram representation of an 8-point DFT function. 
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Figure IB is a block diagram representation of an 8-point DFT obtained by combining two 
4-point DFTs. 

Figure 1C is a block diagram representation of an 8-point DFT obtained by combining four 
2-point DFTs. 

Figure 2A1 is a schematic representation of a DIT radix-2 butterfly computation. 

Figure 2A2 is a schematic representation of a DIF radix-2 butterfly computation. 

Figure 2B is a signal flow graph illustrating an FFT calculation using 12 radix-2 butterflies 
to compute an 8-point FFT in three stages. 

Figure 2C is a signal flow graph illustrating an FFT calculation using 32 radix-2 butterflies 
to compute a 16-point FFT in 4 stages. 

Figure 2D is a signal flow graph illustrating an FFT calculation using 8 radix-4 butterflies 
to compute a 16-point FFT in two stages. 

Figure 3A is a schematic diagram of a DIF butterfly-processing element for an FFT 
algorithm for use with the present invention. 

Figure 3B is a shorthand representation of a DIF butterfly-processing element for use with 
the present invention. 



WO 01/69424 



-9- 



PCT/US01/07602 



Figure 4A is a schematic diagram of a DIT butterfly-processing element for an FFT 
algorithm for use with the present invention. 

Figure 4B is a shorthand representation of a DIT butterfly-processing element for use with 
the present invention. 

Figure 5A is a schematic representation of a DIF butterfly for use with the present 
invention. 

Figure 5B is a schematic representation of a DIT butterfly for use with the present 
invention. 

Figure 6A is a schematic representation of a radix-8 DIT butterfly-processing element for 
use with the present invention. 

Figure 6B is a schematic representation of a radix-16 DEF butterfly-processing element for 
use with the present invention. 

Figure 7 illustrates the design methodology for dividing sets of input data into subsets for 
use in conjunction with the present invention. 

Figure 8 is a block diagram of a parallel FFT multiprocessor in accordance with the present 
invention. 
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Figure 9 is a block diagram of a generalized alternate embodiment of a parallel FFT 
multiprocessor in accordance with the present invention. 

Figure 10 is a block diagram of a 16-point radix-2 embodiment of a parallel FFT 
multiprocessor in accordance with the present invention using four 4-point parallel 
processors. 

Figure 11 is a block diagram of a 16-point radix-2 embodiment of a parallel FFT 
multiprocessor in accordance with the present invention using two 8-point parallel 
processors. 

Detailed description 

The definition of the DFT is shown in equation (1), x (n) is the input sequence, X (k) is the 
output sequence, N is the transform length and w N is the N 1 * 1 root of unity (w N = e" j2n/N ). 
Both x (ri ) and X (k ) are complex valued number sequences. 

n=N-l 

x (ic)= 1 x (n) wgf .k€[0,N-l] (1). 

n=0 

From equation (1) it can be seen that the computational complexity of the DFT increases as 
the square of the transform length, and thus, becomes expensive for large N. A class of 
FFT methods, known as fast algorithms for DFT computation, is based on the divide-and- 
conquer approach. The principle of the divide-and-conquer approach method is that a large 
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problem is divided into smaller sub-problems that are easier to solve. In the FFT case, the 
division into sub-problems means that the input data x n are divided in subsets on which the 
DFT is computed. Then the DFT of the initial data is reconstructed from these intermediate 
results. If this strategy is applied recursively to the intermediate DFTs, an FFT algorithm is 
obtained. Some of these methods are the Common Factor Algorithm, the Split-Radix 
Algorithm, the Prime Factor Algorithm and the Winograd Fourier Transform Algorithm 
(WFTA). Of particular interest is the in place (or butterfly) computation. 

The main objective of many FFT algorithms is a reduction in computation and particularly 
the reduction in the amount of multiplication. 

PARALLEL PROCESSING: Dividing the input data into subsets 

The communication burden is best addressed at the architecture level where the 
representation of the DFT in terms of its partial DFTs is restructured mathematically. The 
problem resides in finding a mathematical model of the combination phase. The problem is 
to find a way that X = T N x (a vector with n components) can be recovered from r vectors 
that are r times shorter. The first step is to divide x itself. The vector S = (x 0 , xi, .., x n .i) is 
split into r shorter pieces and the subdivision continues in till we obtain subsets in each of 
which we could find a single vector of length r. For the first subdivision, the original set So 
(level or layer 0) is subdivided into r subsets: 

Sio = x ra , Si i = X(m+i), and S| (r _D = X( m + (r-i)) (4). 
The subsets of level or layer 1 could be expressed as: 
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for jj = 0, 1,.., r-l and n = 0, 1, ... ((N/r) - 1). Each piece of these subsets are divided into 
r subsets of length (N/r 2 - 1), therefore, each of the subsets S| 0 , Sn, , Si (r -i) is subdivided 
into r subsets yielding to: 
the subset of Sio will include: 

S200 = X r (rn)> ^210 = Xr <m + 1)> • • • •> $2(1-1)0 = x r(m +(r- 1)) (6), 

Sn would have on level (layer) 2 

S2OI = X(r(rn)+1). ^211 = X(r(nH-l >+!)>•., S 2 ( r -l)l = X(r(m+(r- l)+l) (7), 

and finally S Ur .i) on level (layer) 2 will contain: 

S20(r- 1) = X(r(rn)+(r- |)), S2l( r - 1) = X( r ( m+ l) + < r -l)), ,S2( r - l)(r- 1) = *(r(ni + (r- 1)) +(r- 1)) (8). 

The second subdivision could be represented by: 

S 2j,j 2 =X (r(rn + j 2 j + jJ =x (r 2 n + r l j 2 + r ° j, j (9) ' 

for j, = j 2 = 0, 1,.., r -1, and n = 0, 1, .., ((N/r 2 ) - 1). 

In general, the subsets of the i lh level (layer) could be obtained by: 

S ijl.-ji = X (r i n + r< i - , )j i+ r i - 2 j i . 1+ .... + r°j 1 ) (10) ' 
for i = 1,.., p, j; = j M = ...= 0, 1,.., r-l, p = (log, (N)/lo & ) - 1, and n = 0, 1, ... ((N/r* ) - 1). 

In each subdivided subsets, the coefficients just go alternately into r vectors of length r, x<j, 
Xi, and x (r _ X) known as the radix r core calculation (or kernel computation). Then from 
those vectors we form 

Xo = T r x 0 ,X,=T r x,, .....andX^.i) =T r x (r _,) (11). 



WO 01/69424 



-13- 



PCT/US01/07602 



Those are the r multiplications by the smaller matrix T r . As in the replacement of T N by T r , 
the work has been cut by a factor of r. The central problem is to recover X from the (N/r) 
points vectors Xo, Xi, ... , X (r _ i> and Cooley and Tukey noticed how it could be done for 
radix-2 by introducing the concept of the butterfly computation. Mathematically speaking, 
the successive subdivisions of the initial data act as if the one-dimensional transform is 
expressed by a true r* - dimensional transform. Figures 7, 8 and 9 illustrate the division of 
the input data set into subsets. 

PARALLEL PROCESSING: Finding the mathematical model for the combination phase 

Although there is a substantial gain obtained by decomposing the input data into radix r 
core calculations, as shown in figures 2C and 2D, a severe problem arises in the 
complexity of the communication load. For higher radices, the communication load will 
increase drastically with the increasing amount of data. However, the combination phases 
(where partial DFTs are combined) have certain regularity in their structures. The 
combination of r vectors in each subdivided subset of length (3 ((3 = N/r 1 ) is obtained by 
feeding words of length r into the input of the PE. The k* word at the i th level (layer) is 
given by: 

X(k,i) = (X(o,k)> Xd.k). .., X(( r -i) ( k)) (2), 

where k = 0,..,{3 - 1, i = 1,.., p and p = (log r (N)/log r ) - 1. 

The combined output is a vector of length r|3, in which the I th output of the k th word at the 
(i - 1 ) lh level ( layer) is stored into the memory address location x/ according to: 

X /= /r3 + k (3), 
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for/ = 0,1, ,r-l. 

The advantage of appropriately breaking the DFT in terms of its partial DFTs is that the 
number of multiplications and the number of stages may be reduced. The number of stages 
often corresponds to the amount of global communication and/or memory accesses in 
implementation, and thus, reduction in the number of stages is beneficial. Minimizing the 
computational complexity may be done at the algorithmic level of the design process, 
where the minimization of the communication burden is done at the architectural level. 

A fast Fourier transform processor implemented using parallel multiprocessing in 
accordance with the present invention is shown in figure 8. An input memory 802 is 
provided to receive the N input data points. There are r parallel processors 807A through 
807B, each responsive to 1/r of the input 802 data points driven by the source selection 
matrix 806, to compute r partial DFTs. Each of the r parallel processors 807A through 
807B are implemented in accordance with equation (17) for DIF or equation (19) for DIT. 
The outputs of the r parallel processors 807 A and 807B are driven to the input of the 
combination phase 817 by the combination reading address generator 816. After the 
combination phase 817, the output data is driven to the output memory 819 by the 
combination writing address generator 818. The combination phase 817 is implemented in 
accordance with equation (37), given below, for DIT or DIF. 
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Each parallel processor 807A, 807B, comprises a respective reading address generator 
808 A, 808B, writing address generator 809A, 809B, input memory 810A, 810B, output 
memory 811 A, 81 IB and an arithmetic processing unit 814A, 814B. 

A coefficient memory 804 is provided to store the coefficients that will be needed during 
the FFT calculations. In particular, a request for a coefficient is made by any one (or all) 
of the r arithmetic processing unit 814A through 814B. The coefficients needed by each of 
the arithmetic processing unit 814A, 814B, and the respective internal butterfly 
computation units are responsive to the same instructions and use the same coefficients 
804.. 

Note that the r parallel processors 806 A through 806B, do not need to communicate with 
each other. Accordingly, there is no communication burden between the r parallel 
processors 806A through 806B. Each or several parallel processors may be implemented 
on a separate chip and several chips could be mounted on a circuit board without the need 
to send or receive data to or from the other parallel processors, chips or circuit boards. 

The parallel multiprocessor architecture of the present invention is further generalized in 
the block diagram of figure 9. An input memory 902 is provided to receive the N input 
data points. A plurality of parallel processors 906A, 906B, 906C, 906D, 906E, 906F is 
provided. 

Each parallel processor 906A, 906B, 906C, 906D, 906E, 906F includes an arithmetic unit, 
input and output memories, reading and writing address generators, and each has the same 
general structure and function as the parallel processor 807A, 807B in the block diagram of 
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figure 8. Each parallel processor 906A, 906B, 906C, 906D, 906E, 906F is responsive to a 
subset of the N input data points. The outputs of the parallel processors 906 A, 906B, 
906C, 906D, 906E and 906F are combined in respective combination phases 907A, 907B, 
907C, 907D. Combination phase 907 is equivalent to the combination phase 817 in figure 
8 with different parameters. The output of the combination phase 907 A, 907B is stored in 
the input memories of the chips 908A, 908B, 908C, 908D, 908E, 908F (where each chip is 
an electronic layout forming r parallel processors). The data from the chips 908A, 908B, 
908C, 908D, 908E and 908F is combined in combination phases 909A, 909B, 909C, 909D, 
909E, 909F and stored in the input memories of the circuit boards 912A, 912B, 912C, 
912D, 912E, 912F. A circuit board is an electronic device on which r parallel chips are 
mounted. The data of the circuit boards 912A, 912B, 912C, 912D, 912E and 912F are 
combined in combination phases 913A, 913B and stored in the input memories 914A, 
914B of the sub frames 904A, 904B (where a sub-frame is an electronic component on 
which r parallel circuit boards are mounted). The data of the sub frames 904 A and 904B, is 
combined in combination phase 916 and stored in the output memory of main-frame 918. 
It is noted that the further division of the subset of input data need not be divided in the 
same ratio as the original division. 

A specific embodiment of a 16- point radix-2 DFT utilizing the parallel processing 
architecture of the present invention is shown in figure 10. Four separate and independent 
parallel processors 100A, 100B. 100C, 100D are responsive to the input data points. 
Parallel processor 100A comprises adders 110, 112, 114, 116, 124, 126, 128, 130, and 
multipliers 118, 120. Adders 110, 114 and multiplier 118 form aradix-2 butterfly 122. 
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The outputs of the four parallel processors are combined using four radix-4 butterflies 102, 
104, 106, 108 to form the output data points. 

An alternate embodiment of a 16- point radix-2 DFT utilizing the parallel processing 
architecture of the present invention is shown in figure 11. Two separate and independent 
parallel processors 150A, 150B are responsive to the input data points. Parallel processor 
150A comprises a plurality of adders and multipliers. In particular, adders 152, 154 and 
multiplier 156 form a radix-2 butterfly. Similarly, adders 162, 164 and multiplier 166 form 
a radix-2 butterfly. The outputs of the two parallel processors 150A, 150B are combined 
using eight radix-2 butterflies formed by adders 170 and multipliers 172 to form the output 
data points. 

THE JABER PRODUCT ( * ): Representing the DFT in terms of its partial DFTs. 

For a given r x r square matrix T r and for a given set of r column vectors of length 
X 0(fc) ,X 1(k) , ,X (r . 1)(k) , the Jaber product is defined as the operator * as the following 

operation: 
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T r * 



" u (k) 
Xl (k) 



X 



(r-D(k) 



X 0 (0 ) T 00 +X 1 (0 ) T 01 


+ + X (r-l) (0) T 0(r-l) 


X 0 (1) T 00 + X 1 ( |) T 01 


+ + X (r-l) (1) T 0(r-l) 


X 0(^-i) T 00 + 


+ X (r-l)(^-i) T 0(r-D 


X 0 (0) T 10 + X l ( o) T H 


+ + X (r-l) (0) T l(r-l) 




+ X (r-l)/-fl n^Kr-l) 


X 0 (0) T (r-l)0 + 


+ X (r-l)( 0 ) T (r-l)(r-l) 


X 0 ( /3-i) T (r-l)0 + 


••• +X (r-l)(/3-l) T (r-D(r-l) 



which is equal to: 



(12), 





X °(k) 




X 0 (k) T 00 + - + X (r-l) (k) T 0(r-l) 




T r * 


Xl (k) 




X 0 (k) T 10 + •- + X (r-l) (k) Ti(r-l) 


fork = 0,1,.. .,0-1. (13), 




_ X (r-D(k) _ 




_ X 0 (k) T (r-l)0 + •• + X (r-l) (k) T (r-l)(r-l) 





For simplicity the definition of Jaber product expressed in above two equations 
could be represented as: 



r-1 

2 X m (k) T (/,m) 
m=0 



fork = 0,l,...,j3-landZ=0,l,... M r-l (14). 



Mathematical representation of the DFT in terms of its partial DFTs, 
The operation of a radix r butterfly PE is defined by: 

X = B r xx (15), 

where x = [x 0 , x 1? , x r _i] T is the input vector and X= [ X 0 , Xi, . . 

output vector, B r is the rxr butterfly matrix, which can be expressed as: 



X r _! ]' is the 
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B r = W (r,k,i) xT r = [?r DIF(/,mK,j) ] U 



with 



d , _ w ((/mN/r+N(k/r';/r';) N m) 

B rDEF(/,mJl fkii) " w I 1 'J. 



for the decimation in frequency process, 
and 



in which 



B r = T r xW (r k i) = trDlT(/.m)ft,i) ] ^ 

IKJ 

(19), 



B r DIT(/. m ) (k(i) - W 



DITC/.m)^ 

((/mN/r+N (k/r^^^mr^^N 



for the decimation in time process. ((x» N denotes x modulo N, N (k / r 1 ) is defined as the 
integer part of the division of k and r 1 and W (r , ki i) represents the twiddle factor matrix 
expressed by: 



W(r.k,i) ~ 



W 



(0,k.i) 

0 w 



0 

0 - w 



0 
0 

«r-l),k,i) 



""[ W ('.m) (k>i) ] 



(20), 



((N(4> /ri »N 

with w (/tm)(ki) =w ' (21), 



for / = m, and 0 elsewhere, for the decimation in frequency process and 



w,. m , (k , i) =w ((S(k ' r(B ' i))m ' <n - i,))N (22), 



for I = m, and 0 elsewhere, for the decimation in time process. 

T r is an rxr matrix representing the adder-tree in the butterfly, where 
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w° 




- W° 


w° 


w N/r 


w 2N/r 


_ w (r-DN/r 


w° 


w 2N/r 


w 4N/r 


_ w 2(r-l)N/r 


w° 


w (r-l)N/r 




_ w (r-l) 2 N/r 



= tr ( / >m) ] (23), 



«/m-)) N 

where T (; m) = w r (24), 

for/ = m = 0,..., r- 1. 

Based on the proposition in the previous section, equation 1 may be rewritten as: 



N-l 

X (k) = X x (n) W N = T r 
n=0 



(N/r) -1 

X x (rn) W N 
n=0 
(N/r)-l 

x (rn+l) W N 

n=0 



(N/r)-l 

Z (rn+r-l)k 
x (rn+r-l) w N 

n=0 



(25), 



for k = 0, 1, , (N/r) - 1, and n = 0, 1 N - 1. 

Since, 

,rnk ,.,ok 



W N = W N/r 



equation (25) becomes: 

~(N/r)-l 



X(k)=T r * 



X X (m) W N/r 
n=0 

(N/r)-l 

X X (rn+1) w N/r W N 
n=0 



(N/r)-l 

x (rn+(r-l)) w N/r W N 

n=0 



(26), 



= T r * 



'(N/r) -I 

S X (rn) w N/r 
n=0 
(N/r)-l 

W N E x (rn+l) w N/r 
n=0 



(N/r)-l 

W N Zj X (m+(r-l)) w N/r 

n=0 



(27), 



WO 01/69424 



-21- 



PCT/US01/07602 



which for simplicity may be expressed as: 



x (k )=T r ; 



w 



(N/r)-l 

jik V y «/ nk 
N Zj x (rn+j,) w N/r 

n=0 



(28), 



where for simplification in notation the column vector in equation 27 is set equal to: 
(N/r)-l 



w 



X x (rn) w N/r 
n=0 

k (N/rH nk 
W N S x (rn+l) w N/r 

n=0 
(N/r)-l 

(r-Dk V" w nk 
N Zj x (rn+(r-l)) w N/r 

n=0 



(N/r)-l 



.nk 



W N £ x (rn+j,) w N/r 
n=0 



forji = 0,..,(r-l) (29). 



N-l 

^ X oo = S x (o) w n ^presents the DFT of a sequence of length N, 

n=0 



(N/r)-l 



x pOflo - E 



nk 



(N/r)-l 

^pOoo = A rn w N/r> X pl (k) = X X 

n=0 n=0 



(N/r)~l 

rn+l w N/f ,X p(r-l)(k) = X X m+(r-l) 



w N/r 



n=0 



represent the DFT of sequences of length N/r. As a result, equation (28) can be expressed 



as: 



X (k) =T r ? 



w°X 



P0(k) 



W N X pl M 



w (r-l)k x 

W N A p(r-l) (k) 



= T r * [wJj k X pjl(k) ] (30), 



for j, =0, i t .... t r- L 

Finally equation (30) can be developed according to the definition of Jaber product to yield 
the following expression: 
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X (k) - 



w ° X pO M +w 0+k X pl(k)+ . 
w °XpO (k) -w— X pl(k) 



. + w 



0+(r-l)k, 



L P(r-D{k) 

(r-l)N/r+(r-l)k x 
+ ...+ W A P(r-D (k) 



w ° X po w + - 



. + w 



(r-l) z N/r+(r-l)k, 



k P(r-l) (k) 



(31). 



fork = 0,l,....,B- 1. 
Due to the fact that: 



X (k) =T r S 



0 



0 - 



0 0 - w 



0 
0 

(r-Uk 
N 



^00 



X 



= fr r xW (k) ]?k pjl J (32), 



where W(k, = [w (Zijl) ] (33) 

in which w (ZJ]) =w ((jlk))N (34), 

for / = ji and 0 elsewhere, equation 30 in case where X p j represnts the partial DIF or 

(N/r)-l 

DIT DFT of the sequence ^ x (m+j, ) w N/r ma Y be expressed by: 

n=0 

X(k) = [S r (Combination DIT) ]* ^pji (k) ] ^ ^ 

for the decimation in time process, and 

x -fa 1*1x1 (36) ' 

A (k) - Pr (Combination DIF) J* t^pji^ J 



(N/r)-l 



where X p j represnts the partial DIT or DDF DFT of the sequence x (m+ j, ) w 



nk 

+j,) w N/r 



n=0 



Further simplification results when the Z th transform of the k m word is represented by 



th 



r-1 



- 1 X rt , B w N (37) 

Ji=0 
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for/ = 0, 1,..., r- 1, k = 0, 1, ,0-1, where ((x)) N denotes x modulo N. 

The Z 01 transform is stored at the address memory location given by equation (4). 

PARALLEL MULTIPROCESSING ALGORITHM - Figure 8 

The conceptual key of parallel processing is in the reformulation of the Discrete Fourier 
Transform in terms of its partials DFTs that can be combined in a single combination 
phase. The present reformulation enables the design of a multiprocessor parallel structure 
that utilizes r-processing elements in parallel to implement each of the partial DFT 
computations. For a multiprocessor environment, the parallel processor implementation 
results in a decrease in time-delay for the complete EFT of O(r). The time delay could be 
decreased by a factor of 0(r 2 ) by implementing the butterflies with BPEs as described in 
the appendix hereto, which butterfly computing elements utilizes r complex multipliers in 
parallel to implement each of the butterfly computations. 

A second aspect of the combined parallel multiprocessor and BPE in the appendix is that 
that the combination is especially adapted for parallel multiprocessing environments. If 
each butterfly is implemented using the disclosed PBE, each of the r parallel processors 
would always be executing the same instruction simultaneously, which is very desirable 
for implementation. 

The Bit Reversed Input Ordered Output Parallel Multiprocessing Algorithm 
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The data are fed into the memory blocks of each processing element on level (layer) 1 by 
the address generator given in equation (21) Sjj ( =x ( r i n + j )• f° r J' = ®> *•••» r ~* anc * n = 

0,l,..,(CN/r)-l). 

By denoting X p0(k) ,X pl(k) , X p(r . 1)(k) as the partial DIF or DIT DFTs of 

(N/r)-l (N/r)-l (N/r)-l 

£ x rn w N/r- S x rn+l w N/r — and £ x rn+(r-l) w N/r respectively, words 
n=0 n=0 n=0 

of length r are fed into the input of the parallel multiprocessing processing element (MPE). 
The writing address generator w r stores the 1 th MPE's output X/ of the k" 1 word, 

generated by equation (37) above into the input memory of level (layer) 0 at the address 
location: 

W r = l(N/r) + k (39), 

fork = 0, l i ...,(N/i 2 )-l. 

MULTISTAGE PARALLEL MULTIPROCESSING - figure 9 

Further decomposition of the data in terms of its partial DFTs leads to a multistage parallel 
multiprocessing as shown in figure 9. This type of implementation is feasible for 
applications in which a huge amount of data needs to be processed and where rapidity, 
precision and better resolution are required. 
By adopting the following notation: 
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X (k) - T r 0 * 



(40), 



T r * 



-25- 



r j k (N ' r ^ H 
W N L X (r 0 (r rD n + j r( 

n=0 



nk 

r (0 )) W N/r 0 r rD 



N Vd 



w 



N 



k (N/ror (fQ .,))-l 

X X (r 0 (r (lD . 0 n + )+(r 0 - 1 )) w N/r 0 r (r . o 



nk 



n=0 



which may be represented for simplification in notation by: 



wi? k T, * 



w 



" r J0 

N 



(N/r 0 r rj0 )-l 



A «r ri n + j, )+j ID ) w N/r 0 r 
n=0 " r » * 



nk 



(41), 



wjj» k T r * 



w 



JO 



JO 



(N/r 0 r rj0 )-l 



Z, A (r(r rjo n + j rr )+J„ ) w N/r 0 r rj(j 
n=0 J J 0 



.nk 



(42), 



for j rn =0,1, ,r 0 -l,and j r =0,1, ,r (r } -l 

U r- JO 



JO 



In general the radix r FFT expressed in terms of its i* partial FFTs will be as follow: 



X (k) =T r ? 



w 



j 2 k 
N/r 



...T r * 



w 



(N/r')-l 

S 

n=0 

for j! = = =0,1, ,r-landi = l, ,p 



.nk 



S X (r i n + r^j i+ ,.. + r°j 1 ) W N/r i 



(43), 



The conceptual key of multistage parallel processing is the formulation of the radix-r FFT 
as composed of partial FFTs with identical structures and a systematic means of accessing 
the same multipliers coefficients. The latter architecture enables the design of a multi-layer 
processing element (MPE), in which each layer utilizes an r-processing element (PE) in 
parallel to implement each partial DFT computation. In other words, the multi-layer 
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process acts as if the one-dimensional transform is expressed by a true r 1 - dimensional 
transform, which is feasible in neural networks computing. 

The Radix r Multistage Parallel Multiprocessing Algorithm. 



X( K ), elements of layer 0 are obtained by combining the r - 1 subsets of layer one. Each 
subset of layer one is obtained by combining r - 1 subsets of layer 2 (see tree model of 
figure 7). For the layer i, which contains r l PE implemented in parallel with r complex 
multipliers implemented in parallel as shown in figure 9, (N/r 1 ) elements of data are fed 
into the memory blocks of each processing element by the address generator given in 
equation (10) = x^^,,,,. +f(W)j ^ + ^ ) , 

forjj = = ji = 0, 1 ...,r-l. 

At this level, r" partials FFTs are performed in parallel in absence of any 
communication between these parallel processing elements (i.e., are independent of each 
other) and where they will execute the same instruction simultaneously. 

By denoting X:: = ,X iin = ,X iift = the DIT or DIF partials 

FFTs of the subsets on level i, therefore the subsets on level (layer) i - 1 are obtained by: 

,j»k ,, r | _ d._.r^jik 



X < i - 1 )jl OO-Dflc) " Tr * 



W N/r( i " 1 > Xijo ji W 



JO Ji(k) 



(44), 



X G-l)j'l. 



•J'a-0(k) Tf * 



W N 



W N 



0 

(i-D k 



w 



0 
0 

N 





r x ui- 


• ji(0) 








X 


x ij,.. 


■■■*(!) 




x u,... 


j'(r-l) J 



(45), 
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X (i-DJl j(i-l) (k) " k X W (k,i)]*[ X iji j i(k) _ 



where W (k>0 =[w (/ 



(47), 



(46), 



in which w (/ ^ = w ( ^ r< " ,)jik))N (48), 

for I = jj and 0 elsewhere, 

and finally X (i . 1)jl J( . ^ = [B r Combination ^[^j, j, (k) ] (49), 

forj 0 = 0,j, = =ji = 0, 1, ,r- 1 i = p, p-1, ,1 and/ = 0, 1, ....,r- 1. 

Further simplification could results when the 1 th transform of the k ,h word at the i th stage is 
represented by: 

X (i-Dj. " L x ij, w n 

Ji =o 

and the l ,h output of the word at the i* stage will be stored into the memory block of its 
specific PE of the layer at the address location given by: 

X a-Di. J(i-D ( , k .i) =/( 7 )+k (51) ' 

fork = 0, 1, ...,(N/r (,+1) )-l. 

The subdivision of the input data continues till a radix r core calculation is obtained, which 
will generate a radix r FFT algorithm. 



APPENDIX 



Preferred radix-r butterfly -processing element (BPE) 
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Any radix-r butterfly may be used with the present parallel multiprocessing architecture. 
However the preferred radix-r butterfly is assembled from a plurality of butterfly 
processing elements (BPEs). Each BPE, or engine, has a regular structure that can be 
utilized in an array of butterfly-processing elements (BPEs) each having substantially 
identical structure. 

Mathematical basis of the butterfly implementation 



A mathematical term that is a function of r input points and provides a single output point 
is the basis for the design of the present BPE. To provide the insight forming the basis of 
the present BPE, the basic DFT equation is factored to group the variables used in 
multiplication (and simultaneously accessed from memory) into one matrix. In particular, 
starting from the basic DFT equations, the adder matrix is combined with the twiddle 
matrix to form a single phase of calculation. By grouping all the multiply calculations into 
one calculation phase and all the addition calculations into the remaining calculation 
phases, the total number of calculations is reduced and the degree of parallelism is 
increased. 

For a radix-r DIF butterfly, r identical BPEs are arranged in parallel. Each of the r identical 
BPEs are substantially identical to each other and are operated in parallel using the same 
instructions and accessing the necessary set of multiplier coefficients from memory at the 
same time. The outputs of the r identical BPEs form the DFT's r output points. 
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For a radix-r DIT butterfly, (r-1) identical BPEs are arranged in parallel. Each of the (r-1) 
identical BPEs is substantially identical to the others and operates in parallel using the 
same instructions and accessing the necessary set of multiplier constants from memory at 
the same time. The outputs of the (r-1) identical BPEs form the DFT as (r-1) of the r output 
points of the butterfly. The remaining output point (Xo) of the DFT is formed as the sum of 
ther input points. 

Trivial multiplications encountered during the execution of particular butterflies may be 
avoided by simple checks on the coefficient addresses. Avoiding trivial multiplications 
reduces the computational load of particular butterflies. 

An FFT implementation is composed of a plurality of radix-r butterflies with identical 
BPEs and a systematic addressing' scheme for accessing the corresponding multiplier 
coefficients. Each radix-r butterfly utilizes the basic computing unit (BPE), with r (or r-1) 
complex multipliers in parallel to implement each of the butterfly computations* output. 
There is a simple mapping relationship from the three indices (FFT stage, radix-r butterfly, 
butterfly-processing element) to the addresses of the needed multiplier coefficients. The 
simple mapping from the three indices to the addresses of the necessary multiplier 
coefficients accommodates the complexity of higher order radix and mixed radix butterfly 
implementations of the DFT. In a multiprocessor environment, much of the calculations 
are performed in parallel to greatly increase processing speed. 

By using the disclosed BPE in the implementation of the radix-r butterfly, an FFT 
implementation is achieved with a reduced number of calculations and a reduced number 
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of stages of calculations. In addition, the amount of parallelism, both within the butterfly 
processing element (BPE) calculation phases and within the overall FFT algorithm 
butterfly stages permits the use of parallel processing to increase overall FFT calculation 
speed. 

The basic operation of a radix-r butterfly PE is the so-called butterfly in which r inputs are 
combined to give the r outputs via the operation: 

X = B r x x (2A), 
where x = [x (0 ), x ( i), x (r .i)] T is the input vector and X= [ X (0 ), X ( ]), X (r _i) ] T is the 
output vector. 

A special operator, B r is the rxr butterfly matrix, which can be expressed as 

B r =W;xT r (3A) 
for the decimation in frequency process, and 

B r =T r xW^ (4A) 

for the decimation in time process. 

=diag(l, w^j , w^f, ,w ( n' I)p ) represents the twiddle factor and T r is an rxr 

matrix representing the adder-tree in the butterfly, where 

= Er ( /,nJ (5A), 

((/m-)) N 

where T (/ m) = w ' (6A), 



o 


w° 




- w° 


o 


w N/r 


w 2N/r 


_ w (r-l)N/r 


o 


w 2N/r 


w 4N/r 


_ w 2(r-l)N/r 


o 


w (r-l)N/r 




._ w (r-l) 2 N/r 
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1 = m = 0,..., r- 1 and ((x)) N = x modulo N. 



A higher radix butterfly will decrease the number of multiplications and the number of 
butterfly stages, which is the motivation for to implement a higher radix butterfly. Since 
the higher radix reduces the communication load by reducing the number of stages, the 
remaining problem is reduction of the computational burden. Factoring the adder matrix Tr 
(radix 4 case) yielding to the following expression reduces the computational load 



T 4 = 



11 11" 

1 -j -1 j 

1-1 1-1 

1 j -1 -j 



10 1 0' 

0 1 0 -j 
10-10 

0 1 0 j 



X 



10 10 
10-10 
0 10 1 
0 10-1 



(7A). 



As a result, there is a reduction of 4 complex additions for the radix-4 butterfly. For the 
radix-8 butterfly, there is a reduction of 32 complex additions, 4 real additions and 4 real 
multiplications. The reduction of calculations is due to the detection of trivial 
multiplications (i.e., by 1 or -1) and by incorporating (absorbing) the multiplications by j 
or -j in the addition operation by switching the real and imaginary parts of the data. As 
indicated, factoring the adder matrix reduces the number of required computations. 
Although the computation burden is reduced, the complexity of the implementation of such 
butterflies will increase for higher radices. Also, for higher radices, the number of non- 
trivial multiplications increases. 



An important observation is that both the adder matrix T r and the twiddle matrix , 
contain twiddle factors. To control the variation of the twiddle factor during the calculation 
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of a complete FFT, the twiddle factors and the adder matrix are incorporated in a single- 
phase of calculation. 

The special operator B r is defined according to equation (3A), where B r is the product of 
the twiddle factor matrix Wf, and the adder matrix T,. 



So, by defining W (r> kii) the set of the twiddle factor matrices w£, as: 



W, 



(r.k.i) 



w (0,k,i) 0 



0 w 



(l,k,i) 



0 



- w 



0 
0 

((r-l),k,i) 



= L w (/,mKi) 1 



m*k,i) J (8A), 



in which, 



m^tr 1 )) N 



W ('.m) (k , i) - W 



for 1 = m, and 0 elsewhere (9 A), 



therefore, the radix - r computation B r oiF will be expressed: 



to _ w xT _ L 1 (10A), 

u rDIF - w (r,k,i) xl r - PrDIF(/,m| k>i) J 



with B r 



_ w ((/mN/r+S(kfr')lr i )) N 



(HA), 



'rDIF(/,m)( kt i) 

1 = m = 0,. . r - 1, i = 0,1... n - 1, k = 0,1, (N/r) - 1, ((x)) N denotes x modulo N and N (k / 

r 1 ) is defined as the integer part of the division of k and r 1 . 

As a result, the operation of a radix-r PE for the DIF FFT will formulated by: 

thecolonvectorX (r|k)i) =B rDIF xx= [xf/j^J (12A) ' 



r-l 



where the 1 th output X (;) = X x (m )W «"» N " tS (i 3A ), 



m = 0 
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Using the same reasoning as for the radix-r DTP FFT above, the operation of a radix-r DIT 
FFT is derived below. In accordance with equation (4A), B r is the product of the adder 
matrix T r and the twiddle factor matrix , which is equal to: 



B rDIT =T r xW (r,k,i) = l?r 



DIT(Z,m; fk>i) 



j (14A), 



, . , «'mN*+Sck* 0, " i) )nir<"- i) )) N 

in which B rDIT(/ m)(M =w (15A), 



and W (r,k,i) = 



w 



(0,k,i) 
0 w 



(l,k,i) 



0 - w 



0 
0 

((r-D.k,i) 



= [ w (^m; fk>i) ] 



(16A), 



where w (< m)(k = w «So=/r (a -'>)mr( n -'>)) N for j _ ^ Q dsewhere (17A)> 
and n = (log N/log r) - 1. 

As a result, the operation of a radix - r PE for the DIT FFT will be: 

r I (18A) 

thecolonvectorX (rik>i) =B rD1T xx= K/JKiJ 

r-l (((mN/r+fi(k/i (n "' ) )mr (n ' i) )J N 
where the 1 th output X ( t ^ fc j} = J) x (m) w 

m=0 



(19A). 



Radix-r butterfly for an FFT algorithm 



The conceptual key to forming an FFT algorithm is building radix-r FFT butterflies 
composed of BPEs having identical structures and a systematic means of accessing the 
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corresponding multiplier coefficients stored in memory. The present approach reduces the 
number of required multipliers and adders. 

Each butterfly-processing element utilizes either r complex multipliers in parallel (for a 
DIF implementation) or (r - 1) complex multipliers in parallel (for a DIT implementation) 
to implement each of the butterfly computations. There is a simple mapping from the three 
indices (FFT stage, butterfly, processing element) to the addresses of the multiplier 
coefficients needed in figures 3A and 4A. 

Figures 3A and 3B illustrate a DIF BPE for an FFT. A plurality of r complex multipliers 
10A, 10B, 10C, 10D and an adder 14 are provided. One input of each multiplier 10A, 10B, 
10C, 10D is coupled to each respective input data point of the input sequence, x. The 
outputs of each multiplier, 10A, 10B, 10C and 10D are coupled to the inputs of adder 14. 
The output of adder 14 is coupled to one of the output data points of the output sequence, 
X. As shown in figure 3B, a generalized butterfly-processing element 12A, 12B, 12C, 12D 
and 16 is represented by the specialized operator, B r dif- 

Figures 4A and 4B illustrate a DIT butterfly-processing element (BPE) for an FFT. A 
plurality of r complex multipliers 20A, 20B, 20C, and an adder 24 are provided. One input 
of each multiplier 20A, 20B, 20C, is coupled to each respective input data point of the 
input sequence, x. The outputs of each multiplier, 20A, 20B and 20C are coupled to the 
inputs of adder 24. The output of adder 24 is coupled to one of the output data points of the 
output sequence, X. As shown in figure 4B, the specialized operator, B f Drr , represents a 
generalized butterfly-processing element 22A, 22B, 22C and 26. 
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Figures 5A and 5B illustrate how each of the butterfly-processing elements (BPEs) 
represented by the specialized operators, B D if and B D it are arranged to form DIF and DIT 
butterflies respectively. In particular, figure 5A illustrates a DIF radix-r butterfly. In figure 
5A, r butterfly-processing elements 30A, 30B, 30C are arranged in parallel. The inputs of 
each of the r BPEs 30A, 30B, 30C are respectively coupled to each respective input data 
point of the input sequence, x. The outputs of each of the r BPEs 30 A, 30B, 30C are 
respectively coupled to each respective output data point of the output sequence, X. Each 
of the butterfly-processing elements (BPEs) 30A, 30B, 30C, which are operated in parallel, 
are substantially identical as shown in figure 4B. 

Figure 5B illustrates a radix-r DIT butterfly. In figure 5B, r-1 butterfly-processing, 
elements 32A, 32B are arranged in parallel. The inputs of each of the r-1 BPEs 32 A, 32B 
are respectively coupled to each respective input data point of the input sequence, x. The 
outputs of each of the r-1 BPEs 32A, 32B are respectively coupled to each respective 
output data point of the output sequence, X, except for X 0 . Each of the BPEs 32A, 32B, 
which are operated in parallel, is substantially identical to the others as shown in figure 3B. 
Furthermore, an adder 34 is provided to generate the lower order output data point, Xo. The 
r input terminals of adder 34 are respectively coupled to each respective input data point of 
the input sequence, x. The output of adder 34 is coupled to the lower order output data 
point, Xo. 

As shown in figures 5 A and 5B, the present radix - r FFT butterfly is composed of parallel 
computing elements. During the FFT calculation, the simultaneous execution of r 
butterflies in parallel on r separate processors is feasible during each FFT stage. Each the r 
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parallel processors would always be executing the same instruction simultaneously, which 
is very desirable for implementation on some of the latest DSP cards. Trivial 
multiplications encountered during the execution of a particular butterflies (specific 
radices) may be avoided by simple checks on the coefficients addresses. Avoiding trivial 
multiplications reduces the computational load of particular butterflies. Trivial 
multiplications could also be controlled by hardware implementation by using a complex 
multiplier bypass for p =0, where the trivial multiplication could be computed by two 
simple multipliers. 

Reduction in computation can be achieved, by controlling the trivial multiplication in each 
word when (( /mN/r + N(k / r i ) /r i )) N or (( /mN/r M + N(k / r (n *°) mr r i )) N are equal to a 
multiple of N/r. As a result, the most significant advantage is that the multiplication by j 
and -j will be absorbed, yielding to a huge reduction in shifting and negating processes, 
and the only trivial multiplication left is the multiplication by 1 or -1. In this case the 
multiplier is by passed whenever the operators fa) satisfies the condition of being a 
multiple N/r. 

Figures 6A and 6B illustrate the regular internal structure of the present radix-r butterfly- 
processing element (BPE). A radix-8 DIT butterfly-processing element illustrated in figure 
6A includes a plurality of 7 (r-l=7, where r = 8) multipliers 48 and adders 49. The 
multipliers 48 correspond to the (r-1) multipliers 20A, 20B, 20C in figure 4A, and the 
adders 49 correspond to adder 24 in figure 4A. To compute each output data point, there 
are 4 phases of computation. In the first phase 40, multipliers 48 multiply each input data 
point by a coefficient retrieved from memory. Note that all of the coefficients required by 
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multi pliers 48 are retrieved within phase 1. In the next three phases 42, 44, 46, the outputs 
of the multipliers 48 are summed in adders 49, in pairs. That is, in the second phase 42, 
four pair of multiplier 48 outputs are added together. In the third phase 44 two pair of 
adder outputs from the second phase 42 are added together. In the fourth phase 46, one pair 
of adder outputs from the third phase 44 are added together to form a single output data 
point. 

A radix-16 DEF butterfly-processing element illustrated in figure 6B includes a plurality of 
16 multipliers 58 and adders 59. The 16 multipliers 58 correspond to the r multipliers 10A, 
10B, 10C, 10D in figure 3A and the adders 59 correspond to adder 24 in figure 3A. To 
compute each output data point, there are 5 phases of computation. In the first phase 50, 
multipliers 58 multiply each input data point by a coefficient retrieved from memory. Note 
that all of the coefficients required by multipliers 58 are retrieved within phase 1. In the 
next four phases 52, 54, 56, 60, the outputs of the multipliers 58 are summed in adders 59, 
in pairs. That is, in the second phase 52, eight pair of multiplier 58 outputs are added 
together. In the third phase 54 four pair of adder outputs from the second phase 52 are 
added together. In the fourth phase 56, two pair of adder outputs from the third phase 54 
are added together. In the fifth phase 60, one pair of adder outputs from the fourth phase 56 
are added together to form a single output data point. 



As can be seen from figures 6A and 6B, multipliers 48, 58 may be implemented with 
parallel computing elements because the multiplications are occurring simultaneously in a 
single phase. Also, due to the regular structure of the DIT and DIF butterfly-processing 
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elements, an increase in radix, is accommodated by an increase in the number of additional 
addition phases. 

Thus, in applications where minimizing the overall FFT processing time is the most 
important factor, such as in real time signal processing applications, a parallel 
implementation of the radix-r butterfly is beneficial. 
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What is claimed is: 

1: A Fourier transform processor responsive to N input data points to compute N data 
output points, said Fourier transform processor comprising: 



a plurality of r parallel processing elements, each of said plurality of r parallel processing 
elements having N/r input data points and N/r output data points, said N/r input data points 
of said plurality of r parallel processors being coupled to respective ones of said N input 
data points, each of said plurality of r parallel processing elements being independent of 
each other and having at least two stages of butterfly computing elements; and 



a combination phase computing element having N output points forming said N data 
output data points of said Fourier transform processor, said combination phase computing 
element having N input data points coupled to respective ones of said N/r output data 
points from said plurality of r processing elements, said combination phase computing 
element including a single stage of butterfly computing elements. 



2. A Fourier transform processor in accordance with claim 1, wherein said single stage of 
butterfly computing elements in said combination phase computing element is represented 
by: 

X* = X X P Jl(k) w n 07) 
ji=0 

forZ = 0, 1,..., r- 1, k = 0, 1, , (3 - 1, where ((x)) N denotes x modulo N. 
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3. A Fourier transform processor in accordance with claim 1, wherein each of said r 
parallel processing elements includes a butterfly processing element for use in a 
decimation in frequency radix-r butterfly responsive to a plurality of input data points, said 
radix-r butterfly providing a plurality of output data points, said butterfly processing 
element comprising: 

a plurality of multipliers, each of said multipliers having a respective first and second input 
terminals and a respective output terminal; and 

an adder having a plurality of input terminals and an output terminal, each of said plurality 
of input terminals of said adder being coupled to a respective one of said output terminals 
of said plurality of multipliers; 

said plurality of first input terminals of said plurality of multipliers each being respectively 
coupled to one of said plurality of input data points, said output terminal of said adder 
being coupled to one of said plurality of output data points 

said plurality of second input terminals of said plurality of multipliers being coupled to a 
set of coefficients derived from the product of an adder matrix T r and a twiddle factor 

matrix Wn . 
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4. A butterfly processing element in accordance with claim 3, wherein said coefficients 
derived from the product of said adder matrix T r and said twiddle factor matrix wfc , are 
substantially given by 



w° 


w° 


W° 


- W° 


w° 


w N/r 


w 2N/r 


_ w (r-DN/r 


w° 


w 2N/r 


w 4N/r 


_ w 2(r-l)N/r 


w° 


w (r-l)N/r 




_ w (r-l) 2 N/r 



= b(l,m)] ( 5A >> 



((/my))N 

where T (im) = w ' (6A), 



and defining W (ri ^ $ the set of the twiddle factor matrices Wft as: 



W, 



(r,k,i) 



w 



(0,k,i) 



0 



0 w (l,k,i) 

0 0 



- w 



0 
0 

((r-D,k,i) 



= [ w a.m; (kii) ] 



(8A), 



in which, 



(<N(-j)fr')) N 

w ( , m)(kj) =w r for 1 = m, and 0 elsewhere (9 A), 



wherein, 



B r DIF = W (r, k, i) xT r ~ [ B r 



] (10A), 



with B rDIF( , im)(k . ) =w 



DIF(/,m^ kii) 

_ H ,(|/mN/r + (t(k* i )/r i )) N (11A) 
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1 = m = 0 r- 1, i = 0,1.., n - 1, k = 0,1... (N/r)- 1, ((x)) N denotes x modulo N andN (k / 

r") is defined as the integer part of the division of k and r', whereby the operation of said 
butterfly processing element for said decimation in frequency radix-r butterfly is: 

the colon vector X (rki) =B rDEF xx = [fy/j k>i) ] (12A) ' 
where the 1 th output X (IJm = £x (m) w «'»n*+Ro*V» n (13a) 



m = 0 



5. A Fourier transform processor in accordance with claim 1, wherein each of said r 
parallel processing elements includes a butterfly processing element for use in a 
decimation in time radix-r butterfly responsive to a plurality of r input data points, said 
radix-r butterfly providing a plurality of output data points, said butterfly processing 
element comprising; 

a plurality of r-l multipliers, each of said r-l multipliers having a respective first and 
second input terminals and a respective output terminal; and 



an adder having a plurality of r input terminals and an output terminal, one of said plurality 
of r input terminals of said adder being coupled to a first input data point of said plurality 
of r input data points, each of the r-l others of said plurality of r input terminals of said 
adder being coupled to a respective one of said output terminals of said plurality of r-l 
multipliers; 
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said r-1 plurality of first input terminals of said r-1 plurality of multipliers each being 
respectively coupled to one of r-1 of said plurality of r input data points, said output 
terminal of said adder being coupled to one of said plurality of output data points; and 

said plurality of second input terminals of said r-1 plurality of multipliers being coupled to 
a set of coefficients derived from the product of an adder matrix T r and a twiddle, factor 
matrix w£, . 



6. An apparatus in accordance with claim 5, wherein said coefficients derived from the 
product of said adder matrix T r and said twiddle factor matrix w£, , are substantially given 



by, 



B rDIT = T r xW (r,k,i) = tr 



(r,k,i) - PrDIT(/,m| k i) 



(14A), 



((/mN/f+fJOc/r^jmrMMN ' 

in which B rDrr( , m)(k . } = w (15A), 



and w (r,k.i) = 



w (0,k,i) 0 



0 w 



0 



(l,k,i) 



0 



- w 



0 
0 



«r-D.k.i) 



= [ w (/.m| k>i) ] 



(16A), 



where w (/ m)(k j} = w ««** (d '°)»» < "- ,) »n for i = m , and 0 elsewhere (17A), 

and n = (log N/log r) - 1, whereby the operation of said butterfly processing element for 
said decimation in time radix - r butterfly is: 

the colon vector X (rk>i) =B rD i T xx = l*</* kfi) ] (18A) ' 
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r - 1 ((/ m N/r+N (k/A n ' {) )m r (n_i) ;) N 

where the 1 th output X ( = X x (m) w (19A). 

m = 0 

7. A Fourier transform processor in accordance with claim 1, wherein each of said r 
parallel processing elements includes a radix-r butterfly responsive to a plurality of input 
data points and having a plurality of calculation phases, said radix-r butterfly further 
providing a plurality of output data points, said radix-r butterfly comprising: 

a first plurality of multipliers; and 

a second plurality of adders, wherein 

said first plurality of multipliers operate substantially simultaneously during a given 
calculation phase of said plurality of calculation phases, and said second plurality of adders 
operate during said plurality of calculation phases other than said given calculation phase. 

8. An apparatus in accordance with claim 7 wherein said given calculation phase is the first 
calculation phase. 

9. An apparatus in accordance with claim 7, wherein said plurality of multipliers are 
coupled to a set of coefficients derived from the product of an adder matrix T r and a 
twiddle factor matrix . 
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10. An apparatus in accordance with claim 9, wherein said coefficients derived from the 
product of said adder matrix T r and said twiddle factor matrix wfc , are substantially given 
by 





w° 


w° 


W° - w° 




w° 


w N/r 


w 2N/r _ w (r-l)N/r 


T r = 


w° 


w 2N/r 


w 4N/r _ w 2(r-l)N/r 




w° 


w (r-l)N/r 


. _ w (r-D 2 N/r 



= fr(/.mj (5A). 



«'m-)) N 



where T WtB) =w ' (6A), 



and defining W( r , k, o the set of the twiddle factor matrices as: 



W, 



(r,k,i) - 



(0,k,i) 
0 w 



(l,k,i) 



0 
0 



0 



0 " w ((r-l),k,i) 



= I W a.ni> k ,i) 1 



Hi) J (8 A), 



in which, 



w (/ m)(k = w r for 1 = m, and 0 elsewhere (9 A), 



wherein, 



B rDIF = W (r,k.i) xT r= [ B rDIF(/,m| kii) ] 



with B rDnS(< m)ai) =w 
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1 = m = 0,..., r - 1, i = 0,1.., n - 1, k = 0,1.., (N/r) - 1, ((x)) N denotes x modulo N and N (k / 
r') is defined as the integer part of the division of k and r', whereby the operation of said 
butterfly processing element for said decimation in frequency radix-r butterfly is: 

thecolon vector X (r k ;) =B rDIF xx = [fy/^ ki) ] ^ 12A ^' 



r-l 



where the 1 th output X (< = £x (m) w 



((/mN/r+N(k/r')/r')) N 



(13A). 



n = 0 



11. An apparatus in accordance with claim 9, wherein said coefficients derived from the 
product of said adder matrix T r and said twiddle factor matrix w£, , are substantially given 
by, 



B rDIT =T r xW (r,k,i) = PrDrT(/,mJf k)i) 



(14A), 



in which B 



rDIT(/,m) (ki i) 



((/mN/rtf!(Wr M )mr M ))N „ , . . 

= w (15A), 



and W (r k i) = 



w 



(0,k,i) 
0 w 



0 



0 

(U,i) 



0 



0 
0 

w ((r-l),k,i) 



= [ w (i,m* k(i) ] 



(16A), 



where w ( , . } = w «*»o* <0 " i) )' nr(n " i) »N for \ = m> and 0 elsewhere (17A), 

and n = (log N/log r) - 1, whereby the operation of said butterfly processing element for 
said decimation in time radix - r butterfly is: 

the colon vector X (rki) =B rDIT xx = ^/^ ki) ] (18A) ' 
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r-1 

where the l lh output X(i) (k i) = £*(m) w 



((!mN/r+N(k/r (n,i) ;mr (n_i) JJ N 



(19A). 



m = 0 



12. A multistage Fourier transform processor responsive to N input data points to compute 
N output data points, said multistage Fourier transform processor comprising: 

a first multiprocessor stage coupled to said N input data points, said first multiprocessor 
stage comprising a first plurality of Fourier transform processors, each of said first 
plurality of Fourier transform processors comprising 

a plurality of rl parallel processing elements, each of said plurality of rl parallel 
processing elements having N/rl input data points and N/rl output data points, said 
N/rl input data points of said plurality of rl parallel processors being coupled to 
respective ones of said N input data points, each of said plurality of rl parallel 
processing elements being independent of each other and having at least two stages 
of butterfly computing elements; and 

a combination phase computing element having N output points, said combination 
phase computing element having N input data points coupled to respective ones of 
said N/rl output data points from said plurality of rl processing elements, said 
combination phase computing element including a single stage of butterfly 
computing elements; and 
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a second multiprocessor stage coupled to said N output data points of said first 
multiprocessor stage, said second multiprocessor stage comprising a second plurality of 
Fourier transform processors, each of said second plurality of Fourier transform processors 
comprising, 

a plurality of r2 parallel processing elements, each of said plurality of r2 parallel 
processing elements having N/r2 input data points and N/r2 output data points, said 
N/r2 input data points of said plurality of r2 parallel processors being coupled to 
respective ones of said N output data points of said first multiprocessor stage, each 
of said plurality of r2 parallel processing elements being independent of each other 
and having at least two stages of butterfly computing elements; and 

a combination phase computing element having N output points forming said N 
data output data points of said Fourier transform processor, said combination phase 
computing element having N input data points coupled to respective ones of said 
N/r2 output data points from said plurality of r2 processing elements, said 
combination phase computing element including a single stage of butterfly 
computing elements. 

13. A multistage Fourier transform processor in accordance with claim 12, further 
including a plurality of (i) multiprocessor stages, wherein the radix r FFT of each said 
plurality of (i) multiprocessor stages is expressed in terms of its i th partial FFTs as follows: 
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X (k) =T r i 



w jj k T r * 



w 



j 2 k 
N/r 



...T r * 



w 



N/r< M > 



(N/r/M 

S X (r i n + r< i - , >j i ^>r°j l ) W 



nk 
N/r 1 



(43), 



n=0 



for jj = = ji = 0 f 1, , r - 1 and i = 1, , p. 
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FIG. 1A 
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FIG. 2A1 FIG. 2A2 
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FIG. 2C 
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FIG. 3A 
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FIG. 6 A 
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FIG. 8 
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FIG. 9 
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